d <- read.csv("sweep21/sweep21.frequencies",header=FALSE)

dat <- d[d[,3]>0.05&d[,1]>0,]

rrr <- unique(dat[,1])
rrrc <- gray(rrr/max(rrr))
clonecol <- cbind(rrr,rrrc)

plot(dat[,2],dat[,3],col="black",ylab="Frequency",xlab="Time (days)",pch=19)
for(i in 1:length(rrr)){
lines(dat[dat[,1]==rrr[i],2],dat[dat[,1]==rrr[i],3],col="black")
}

# Color by clone selective advantage

# Load clones selective advantage, color clones by advantage
clones.fitness <- read.csv("sweep21/sweep21.clones",header=FALSE)
d <- read.csv("sweep21/sweep21.frequencies",header=FALSE)

# Remove duplicates
clones.fitness <- unique(clones.fitness)

cfit <- cbind(rep(0,nrow(clones.fitness)),rep(0,nrow(clones.fitness)))
cfit[clones.fitness[,2]==1,1] <- "black"
cfit[clones.fitness[,2]==2,1] <- "orange"
cfit[clones.fitness[,2]==4,1] <- "red"
cfit[clones.fitness[,2]==8,1] <- "brown"
cfit[clones.fitness[,3]==1,2] <- "black"
cfit[clones.fitness[,3]==2,2] <- "orange"
cfit[clones.fitness[,3]==4,2] <- "red"
cfit[clones.fitness[,3]==8,2] <- "brown"

# get clones that have reached more than 5% frequency over 20 years
# clone_id, time, frequency
dat <- d[d[,3]>0.05&d[,1]>0,]

intersect <- function(x, y) y[match(x, y, nomatch = 0)]

# get the frequency of clones in dat, even if it dropped to 0 at some time points
clones.freq <- d[d[,1] %in% dat[,1],]

# Plot
par(mfrow=c(1,2))
plot(clones.freq[,2],clones.freq[,3],col="black",ylab="Frequency",xlab="Time (days)",main="Survival Advantage",pch=NA)
rrr <- unique(clones.freq[,1])
for(i in 1:length(rrr)){
lines(clones.freq[clones.freq[,1]==rrr[i],2],clones.freq[clones.freq[,1]==rrr[i],3],col=cfit[clones.fitness[,1]==rrr[i],1])
}
plot(clones.freq[,2],clones.freq[,3],col="black",ylab="Frequency",xlab="Time (days)",main="Reproductive Advantage",pch=NA)
rrr <- unique(clones.freq[,1])
for(i in 1:length(rrr)){
lines(clones.freq[clones.freq[,1]==rrr[i],2],clones.freq[clones.freq[,1]==rrr[i],3],col=cfit[clones.fitness[,1]==rrr[i],2])
}



# Load data
clonesn21 <- read.csv("sweep21/sweep21.clones.number.log",header=FALSE)

par(mfrow=c(1,1))
plot(clonesn21[,1],clonesn21[,2],col="black",ylab="Absolute number of clones",xlab="Time (days)",main="10e-06",type="l")




# Plot number of clones over time
# Load data
clonesn6 <- read.csv("sweep3/sweep3.clones.number.log",header=FALSE)
clonesn7 <- read.csv("sweep5/sweep5.clones.number.log",header=FALSE)
clonesn8 <- read.csv("sweep6/sweep6.clones.number.log",header=FALSE)

par(mfrow=c(1,3))
plot(clonesn6[,1],clonesn6[,2],col="black",ylab="Absolute number of clones",xlab="Time (days)",main="10e-06",type="l")
plot(clonesn7[,1],clonesn7[,2],col="black",ylab="Absolute number of clones",xlab="Time (days)",main="10e-07",type="l")
plot(clonesn8[,1],clonesn8[,2],col="black",ylab="Absolute number of clones",xlab="Time (days)",main="10e-08",type="l")

grid6 <- read.csv("sweep6/sweep6.log",header=FALSE)
ggrid6 <- grid6[,2:3]
cnum <- NULL
ctime <- NULL
for(i in 1:40){
	cnum[i] <- nrow(unique(ggrid6[ggrid6[,2]==i*182,1:2]))
	ctime[i] <- i*182
}
par(mfrow=c(1,1))
plot(ctime,cnum,col="black",ylab="Absolute number of clones in grid",xlab="Time (days)",main="10e-08",type="l")








